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Abstract 

We study an interacting ensemble of instantons at finite baryon chemical 
potential. We emphasize the importance of fermionic zero modes and calcu- 
late the fermion induced interaction between instantons at non-zero chemical 
potential. We show that unquenched simulations of the instanton ensemble 
are feasible in two regimes, for sufficiently small and for very large chemi- 
cal potential. At very large chemical potential chiral symmetry is restored 
and the instanton ensemble is dominated by strongly correlated chain-like 
configurations. 
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I. INTRODUCTION 



The behavior of matter under extreme conditions of temperature and/or density is of 
great interest in connection with the physics of heavy ion colhsions and neutron stars. In 
addition to that, understanding the nature of the phase transition between a hadronic gas at 
low temperature and density and the quark gluon plasma at high T and pB is an important 
problem in its own right and expected to shed some light on the vacuum structure of QCD. 

While the structure of QCD at high temperature has been investigated in some detail 
using lattice simulations little is known about matter at high baryon density. This 
is unfortunate, because the possible phase structure of matter at finite density is much 
richer, including phenomena such as pion P| and kaon condensation density isomers 
strange quark matter [^, etc. From the point of view of lattice QCD, the difficulty 
is connected with the sign problem (see P] for a recent review on lattice simulations at 
finite chemical potential). At zero temperature and chemical potential, the weight factor 
det(^ + m)exp(— S") in the euclidean path integral is real and positive. This means that 
the weight function can be interpreted as a probability distribution and one can perform 
simulations based on importance sampling. This remains true also at finite temperature. 
In euclidean space, the heat bath is represented by the (anti)periodic boundary conditions 
imposed on the gauge fields and fermions. At finite density, however, the determinant 
det(-^ + m — /i74) is complex and cannot be interpreted as a probability distribution. 

Furthermore, while it is sometimes argued that neglecting the determinant (the so-called 
quenched approximation) provides a reasonable approximation at finite temperature, the 
same cannot be true at finite chemical potential. If the determinant is neglected, or if 
only the modulus of the determinant is retained, the critical chemical potential for chiral 
symmetry restoration is zero. This phenomenon is easily understood in the case where 
we retain the modulus only [^. This approximation corresponds to a theory where half 
of the quarks couple with the opposite chemical potential. In this case, the lightest state 
that carries baryon number is a baryonic pion, a pion made of an ordinary quark and an 
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anti-quark that couples with the opposite chemical potential. This state is degenerate with 
the pion, so that in the chiral limit it can Bose condense even if the chemical potentials is 
arbitrarily small. In full QCD, on the other hand, the baryon density has to be zero as long 
as /i < — Mb/ 3, where Mb is the mass of the lightest baryon. 

Recently, this phenomenon has been studied in some detail in random matrix models 
of QCD [p|-[To[|. While these models nicely illustrate the shortcomings of the quenched 
approximation they do not really improve our understanding of the phase transition in 
full QCD. The chemical potential is simply added as an external field to the random matrix 
that represents the Dirac operator. Chiral symmetry is restored when the chemical potential 
exceeds a certain critical value, but the dynamics of the phase transition is equivalent to a 
simple mean field model with static quarks [ pJ]) . 

In this work, we want to study QCD at finite chemical potential in the instanton liquid 
model. The instanton model is based on the assumption that the physics of chiral symmetry 
breaking and the dynamics of light hadrons is dominated by just a few low-lying modes in 
the spectrum of the Dirac operator. These modes are linear combinations of the zero modes 
associated with instantons, classical solutions of the euclidean Yang-Mills field equations. 
The instanton model provides a phenomenologically very successful description of the QCD 
vacuum, chiral symmetry breaking and hadronic correlation functions [jl2|. In addition to 



that, the underlying assumptions of the model can be checked in great detail on the lattice, 
see e.g. [pr3| , p!4| . Of particular interest are recent simulations that focus on the correlation 
of low-lying eigenstates of the Dirac operator with topological charges, and on the role that 



these states play in the correlation functions of light hadrons |1T5|- |T^] . There is evidence that 
low-lying states are approximately linear combinations of instanton zero modes, and that 
these states play a dominant role in the correlation functions of light mesons. 

Studying instantons at finite chemical potential is worth while for several reasons. First, 
the model provides a very interesting mechanism for the chiral phase transition at finite 
temperature, based on the formation of topologically neutral instanton-anti-instanton pairs 
near T = Tc [|l^]. The corresponding eigenstates are localized modes, contrary to the 



delocalized states that make up the quark condensate at T = 0. This scenario raises the 
question whether a similar mechanism might be at work at finite baryon density. Second, 
because the model concentrates on just a few modes in the spectrum of the Dirac operator, 
the sign problem should not be as severe as it is in lattice simulations. This should allow 
more realistic simulations to be performed. Finally, the instanton model can be studied in a 
number of approximations that allow for analytical (or almost analytical) solutions. In this 
case, one might be able to handle the analytic continuation required to introduce a chemical 
potential. 

This paper is organized as follows. In Sect. 2 we discuss zero modes, the anomaly and 
the density of instantons at // 7^ 0. In Sect. 3 we calculate matrix elements of the Dirac 
operator between zero mode wave functions. In Sect. 4, we study the spectrum of the Dirac 
operator in an ensemble of instantons at 7^ 0. In Sect. 5, we present simulations of the 
instanton liquid at small and large chemical potential. 



We would like to begin our discussion by constructing the free euclidean coordinate space 
propagators of bosons and fermions at non-zero chemical potential. This means that we are 
looking for solutions of the equations {ida — iSailJ'Y^ix) — 5{x) and {i^ — iix'y^S{x) — 
5{x) subject to the appropriate boundary conditions. Naively, we can obtain solutions by 
multiplying the zero density propagator by exp(//a;4). This means that at finite chemical 
potential the propagation of quarks along the positive time direction (and anti-quarks along 
the negative time direction) is favored over propagation in the opposite direction. 

However, this naive prescription does not provide the correct propagator. We have missed 
the contribution from occupied states below the Fermi surface. The correct solution is most 
easily obtained by Fourier transforming the momentum space propagator. The result for 
the charged Klein- Gordon propagator is 



II. ZERO MODES AND THE ANOMALY AT // / 




(1) 
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Note that the extra contribution from occupied states has cancelled the exponent exp(/iX4). 
From (|1|) we can also obtain the Fermion propagator 

Six) = {%^ + Z/i74) A(x) = (7 ■ fSr{x) + 74^4(X)) • (2) 

where the two functions Sr{x) and 5*4 (x) are given by 

Sr{x) = (2r — 11x^x4^/ cos(/ir) + {x\ + 3r^a;4 + jj.r'^x^^ \ (3) 

5*4(0;) = (^2x4 + ^x"^^ cos(/ir) + 1^X4 — + ^x^x"^^ __(aO ^ ^^■^ 

Asymptotically, the propagator behaves as S{x/i) ~ Ai^/x4 in the temporal direction, as 
compared to ~ I/X4 at yU = 0. The propagator in the positive X4 direction is still enhanced, 
but only by a power of X4 rather than an exponential. This enhancement is related to the 
fact that at finite /i, quark paths dominantly wind around the imaginary time direction. The 
number of loops piercing an X4 = const surface going forward in time, minus the number of 
loops going backwards, is related to the total baryon number. From (^) one easily reproduces 
the perturbative result pb = ^Si{r — > 0,0:4 = 0) = /i^/(37r^). In the next section we will 
also see that the anisotropy of the propagator causes the interaction between instantons to 
be strongly anisotropic as well. 

For spacelike separations the propagator behaves as S{r) ~ /isin(yur)/r^. It decays 
as and oscillates with a wave number given by the Fermi momentum pp = /i. The 
oscillations are due to the presence of a sharp Fermi surface and are analogous to Friedel 
oscillations in a degenerate electron gas . 

Zero modes in the spectrum of the Dirac operator play a central role in understanding 
the U{1)a anomaly and the mechanism for chiral symmetry breaking. While the relation 
to the anomaly is a rigorous consequence of the arguments given below, the connection 
with chiral symmetry breaking is more subtle. The Banks-Casher formula ||22| (qq) = 
—npi^X = 0) relates the density of eigenvalues of the Dirac operator at zero virtuality to the 
quark condensate. In a finite system with volume V these states are not exact zero modes, 
instead they occupy a band around zero with a typical level spacing AA ~ \{<l(l)\/y p3| . 
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Perturbative states, on the other hand, are spaced much farther apart, AA ~ It is 

tempting to identify the extra states with quasi zero modes that arise when the exact zero 
modes of isolated instantons and anti-instantons interact. The fact that this mechanism 
leads to chiral symmetry breaking was checked in simulations of the instanton ensemble 
in QCD Whether quasi zero modes associated with topology dominate the quark 



condensate in QCD is now under active investigation on the lattice [p!6| , p!7tp6| , and initial 
results are very promising. 

The connection of zero modes with the axial anomaly can be seen from the fact that the 
violation of axial charge in an arbitrary gauge potential is given by AQs = 2Nf{nL — tlr), 
where n^^^ is the number of left and right handed zero modes of the Dirac operator. Index 
theorems (or the anomaly equation = N f / {1Qti'^)G ^yG ^y) then imply that zero modes 
have to be associated with topological charges. This requirement is satisfied by the instanton 
solution. The instanton has topological charge Q = +1, and the Dirac operator in the 
instanton field has a left-handed zero mode, ipipQ = and 75'?/^o = —i'o- Analogously, there 
is a right handed zero mode in the field of an anti-instanton. 

It is essential that these states are exact zero modes of the Dirac operator, because every 
non-zero mode automatically occurs in pairs and only unpaired states contribute to AQ^. If 
ipx is an eigenstate of iJp with eigenvalue A 7^ 0, then •y^ipx is an eigenstate with eigenvalue 
—A. This remains true at yU 7^ 0, because the extra term ifi'jA has the same chirality as ip. 
But the anomaly is a short distance effect and should not be affected by a small chemical 
potentia^^. We therefore expect that there should be an exact zero mode in the spectrum of 
the Dirac operator in the instanton field also at /i 7^ 0. 

This state can be found using the following observation ||30|: Naively, a solution of the 
Dirac equation at yU 7^ is given by exp{fiX4^)ipQ{x), where ipo{x) is the /i = zero mode. From 
our discussion of the quark propagator we expect that this is not an acceptable solution. 



^See |27-2S] for a similar discussion at finite temperature. 
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Instead, we have to construct a solution of the /i = equation that cancels the exponential 
growth. 



Such a solution can be found from the general ansatz 31 



where Xi = ^ia is a Dirac and color spinor and 7± = (l±75)/2 for instantons/anti-instantons. 
The gauge potential of the instanton is given by Ay = rj^pdp In n(x) where fj^^^ is the anti- 
self-dual 't Hooft symbol ||32| and n(x) = 1 + p^/x^. For an anti-instanton, we have to 
make the usual replacement r/^^, r]'^^. Substituting the ansatz (^) into the Dirac equation 
leads to the condition that $(a;) has to be a solution of the Laplace equation = 0. 

At /i = 0, a solution that is regular and satisfies the normalization condition is given by 

This suggests that a solution which satisfies the boundary condition at finite chemical 
potential can be found from the /i 7^ Klein-Gordon propagator A{x). If A (a;) satisfies the 
Klein-Gordon equation with a chemical potential, then $(a;) = (47r^)p^ exp(— /ia;4)A(x) is a 
solution of the Laplace equation that behaves as ^{x) ~ exp{—fiX4). It is easy to check that 

M^) = ^/2npexpi^^x,)^W (^J^Pt^^A^\ u'^^^^x' (6) 



indeed solves {ip — ifi^4)ip(j = p0| , |33| . Furthermore, the zero mode solution is well behaved 
at infinity. It decays like the free quark propagator (Q). Some care has to be taken in 
constructing the adjoint solution. The Dirac operator is not hermitean, so left eigenstates 
are not the hermitean conjugate of right eigenstates. In fact, it is easy to see that V'ol/^) ^) = 
[ipo{—fi,x)]* is a solution of ipoiip — ip^yA) = 0. Using this property, we can also check that 
the zero mode is correctly normalized, 

(fxil)l{p,x)il)o{p,x) = l. (7) 

In order to determine the importance of instanton effects at finite baryon density, we 
have to know the dependence of the instanton rate on the baryon chemical potential. At 



large chemical potential, this dependence can be determined from a perturbative calculation. 
The result isQ |3gj3|j3|j3|g. 

n(p, /i) = n{p, 0) exp (-Ar/(p/i)^) , (8) 

where n{p,0) is the /i, T = rate originally calculated by 't Hooft [^]. The suppression 
factor exp(— A^/(p/i)^) is completely analogous to the factor exp{—2{Nc/3 + Nf /6){npT)'^) 
that appears in the rate at finite T In both cases, the physical reason for the suppression 
of instantons is Debye screening of the strong 0(1/ g) color electric fields inside the instanton. 
The rate at T 7^ contains some extra terms that are due to the fact that the instanton 
solution itself is deformed at T 7^ 0. At least in perturbation theory, this effect is absent at 

Since the high temperature suppression factor is due to the scattering of perturbative 
quarks and gluons on the instanton, it was argued that the result should not be trusted 



at temperatures around and below Tc [^. At low temperature the rate is determined by 
the scattering of hadrons, not quarks and gluons, on the instanton. This has been checked 
in lattice calculations of the temperature dependence of the topological susceptibility in 
quenched QCD p9|-p]]. These calculations find the expected exponential suppression above 
Tc, but below Tc the topological susceptibility is essentially independent of T. 

Analogously, we expect that the result (^) is only valid if the chemical potential is 
large. At small baryon density, the rate is determined by the density dependence of the 
effective 2A^j-fermion operator induced by instantons. In the case of small T and zero p 



this dependence can be determined from soft pion theorems |38|. At finite density, only the 



behavior of the quark condensate is known |^ 



{m) = {m)[i-^^PB + ..j, (9) 

where S^tv — 46 MeV is the nN sigma term. This implies (qq) ~ {qq)o{l — 0.3{pb/Pb)) 
where p% = 0.14fm~^ is nuclear matter density. If we estimate the density dependence 



See [34 1 for a discussion of the disagreements among some of the earlier calculations. 
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of the 't Hooft operator using factorization, we get a very rapid dependence, {Ohh) — 
{0'tH)o{^ — ^-^{Pb/ P%)) for Nf = 2. This seems too large, and factorization is not expected 
to be a good approximation at finite density. In order to improve on this estimate we have 
to determine the expectation value of the 't Hooft operator inside the nucleon. We have not 
performed this calculation, but it appears straightforward in principle. In practice, we have 
not used the perturbative result (|^) in our simulations at small chemical potential, but we 
have taken it into account for large /i. 

III. MATRIX ELEMENTS OF THE DIRAC OPERATOR 

In order to study an ensemble of instantons at finite chemical potential, we have to de- 
termine matrix elements of the Dirac operator {ip — ifi'j^) between the zero modes wave 
functions tpj associated with different instantons. At finite temperature, the functional form 
of these matrix elements plays a crucial role in determining the mechanism for the chiral 
phase transition. Near Tc, the overlap matrix elements favor the formation of polarized 
instanton-anti-instanton molecules |4^JT9|] . Our aim here is to identify instanton configura- 
tions that may play a similar role at large chemical potential. 

In the following, we will assume that the gauge potential of a system of instantons 
and anti-instantons can be approximated by the sum of the individual gauge potentials 
A,, = J2i ^i- In this case, the overlap matrix element is given by 

TjA = J d'^x ipiix — z) (ip — ifx'ji) ipAix) = J d'^x^l^x — z) + ifx'ji) ipAi^). (10) 

Here, we have used the equations of motion to reduce the covariant derivative to an ordinary 
one. Also note that because of the chirality of the zero modes, only matrix elements between 
instantons and anti-instantons are non- vanishing. 

The overlap matrix elements Tja_ depend on the distance z^ between the instantons, 
the instanton sizes p/.^ and the relative color orientation U = U\Ua- In the following, we 
will characterize the color orientation using the four- vector Uy = l/(2i)tr([/rj^+^), where 
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rj='=) = (t, =f). At fi,T = 0, Lorentz invariance implies that Tja = i{u ■ p/, p^). 

At p,T 7^ 0, Lorentz invariance is broken and the matrix elements have the more general 



structure (see for a discussion of Tia at T 7^ 0) 

TiA = iUih+i{u- z)f2 (11) 

with fi = fi{\z\, Z4, Pi, Pa) and z = z/\z\. These functions can be calculated by writing the 
zero mode solutions in the form ip1 = 4>u(yX){lu)ijU°'^x^j- At p,T = 0, we have 0jy(x) = 
x^(j){\x\), while in the more general case = 5j,404(r, X4) + 6^ifi4>r{r,Xi). In terms of 0^, 

the overlap matrix elements are given by 



Tja = 2iu^ j d'x {<Ptdu<Pi - <Ptda<Pi + (t>tdAi - P {<PM - Sa4<pUi + <PUi)) , (12) 
where (j)^ = (p^{x — z) and 0^ = 0i/(a;). From (|1^) we can read off the invariant functions 

/l,2 

/i = 2 1 d^x {<Pi [d^<Pi + + dr<pi) + {r^z) ■ r (dr<pi - d^(t>i) (13) 

-p(ct>i<p{-{f^z)-f<pf<Pi)], 
/2 = 2 1 d^x {z ■ (f^z) (f)f (940i + 2<Pi/r + 9,0^ +f-z<Pi (^40^ - dr<Pi) (14) 

- p ■ f + z ■ (r^z) (l)^(pi) } . 

The functions 0r.,4 are easily determined from (^) but the resulting integrals are fairly com- 
plicated and have to be performed numerically. Results are shown in Fig. |l] which also 
displays a simple parametrization of the matrix elements given in App. 0. Similar results 
have been obtained by Rapp [^]. The most important point is that the asymptotics of Tja 



are the same as those of the free quark propagator. In particular, Tja ~ p'^/x^ in the time 
direction and Tja ~ psin(pr)/r^ in the spatial direction. 

In Sect. ^ we saw that at finite baryon density, quark paths tend to wind around the 
time direction. The probability of an instanton configuration is proportional to the fermion 
determinant, which contains the matrix elements Tja- Since these matrix elements have the 
same functional form as the quark propagator, we also expect the instanton ensemble to be 
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dominated by instanton-anti-instanton chains that wind around the time direction. We will 
study this phenomenon in Sec. [V|. 



IV. THE INSTANTON ENSEMBLE AT fi ^ 

The instanton ensemble at finite baryon chemical potential is determined by the partition 
function 



Here, and are the numbers of instantons and anti-instantons, the correspond- 
ing collective coordinates (position, color orientation and size), and n{p) the semi-classical 
size distribution. The bosonic interaction between instantons is denoted by Smt while the 
fermionic interaction is contained in the determinant det{p + m — /i74). 

Following the usual strategy in instanton model we will split the fermion determinant 
into a high and a low momentum part. The low momentum part is related to fermion zero 
modes and will be treated exactly. For the high momentum part we will assume that it 
can be factorized into the contributions of individual instantons. For a single instanton, the 
non-zero mode determinant in Gaussian approximation is included in the size distribution 
n{p). The dynamics of the system is then governed by the determinant in the zero mode 
basis. At finite chemical potential, the Dirac operator has the structure 



where T is the matrix of overlap matrix elements Tja- The block off-diagonal form of the 
Dirac operator is due to the chiral structure of the zero modes. Clearly, for p ^ the Dirac 
operator is not hermitean and its eigenvalues are complex. 

For comparison, the form of the Dirac operator in a basis of /i = zero modes is 




(15) 




T{p) 
^Tt(-/.) 



(16) 



/ 







T(^ = 0) +/iM(/i = 0) 




(17) 



V 



Tt(/i = 0) -/iMt(/i = 0) 
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where Mja are matrix elements of 74. This is the structure of the Dirac operator that is 
assumed in random matrix models of QCD at finite chemical potential (with the additional 
simplification M I) |[. 

We will study the difference between the exact matrix elements (0) and the simplified 
version (^) in Sec. |VI[ But even without a calculation it is clear that the full matrix 
elements can carry significant dynamics whereas in the simplified form (in particular if the 
substitution M — 1 is made) the chemical potential only acts like an external field. The 
situation is similar to what happens at finite temperature. While random matrix models 
provide a good description of the dynamics of zero modes at T = (and Nj not too large), 
the dynamics of the phase transition is very different (see for a random matrix model 



that mimics the phase transition in the instanton model). There is one more important 
difference between the Dirac operators ([161) and (0): For ^ 00 the schematic Dirac 



operator ([TTD becomes anti-hermitean, so all eigenvalues are imaginary and the determinant 
is real. For the exact zero modes, on the other hand, the Dirac operator does not simplify 
in that limit. 

The statistical system described by the partition function (^) is quite complicated so 
that in general one has to rely on numerical simulations to determine the phase structure. 
At finite temperature and zero chemical these simulations are straightforward but at 
finite chemical potential one has to deal with the fact that the determinant is complex. 

At this point, the only practical solution to this problem is to generate an ensemble with 
the absolute magnitude of the determinant 

^ N++N- Nf 

^11= AT \N \ n Win{pi)\exp{-Sint)\X{det{p + mf - ^i-fi) , (18) 

and incorporate the phase information into the observables. This means that expectation 
values are determined from 

(O) = ^ (19) 

where (.)|| is an expectation value determined with the modulus of the determinant and 
is the phase of the determinant. As mentioned in the introduction, taking into account the 
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phase information is absolutely essential. The theory described by the partition function 
with the modulus of the determinant is very different from QCD. In particular, we expect 
chiral symmetry restoration to take place at ~ 

The problem with this method is that if fluctuations in the phase are large then (e**^) 
will typically be close to zero and the expectation value (|T^ is essentially undetermined. In 
other words, if fluctuations are big then the ensemble generated by (p^Sf) has a very small 
overlap with the correct ground state. 

An important simplification occurs in the special case of two colors, Nc = 2. In this case, 
the matrix elements Tja are real (because the color orientation vector is real). Eigenvalues 
not only occur in pairs ±A, but also in complex conjugate pairs A, A*. Not all eigenvalues 
are four-fold degenerate because some of them are purely real or purely imaginary. This 
means that the determinant is real, but not necessarily real and positive. We can construct 
a theory with a real and positive determinant by considering the case Nc = 2 and Nf even. 

Since the determinant is positive, QCD with = 2 and Nf = 2 is easy to simulate even 
if the chemical potential is non-zero. But as in the case of quenched QCD, the physics is 
very different from real QCD. If A^^^ = 2, baryons are diquarks. Furthermore, the N^ = 2 
theory has an additional particle-anti-particle (Pauli-Giirsey) symmetry symmetry 
which implies that the scalar diquark is degenerate with the pion. This means that the 
critical chemical potential is again expected to be yU ~ mj^/2. 



V. SPECTRUM OF THE DIRAC OPERATOR AT FINITE CHEMICAL 

POTENTIAL 

The low virtuality part of the Dirac spectrum is intimately connected with the physics 
of chiral symmetry breaking. For /i = 0, this is expressed by the Banks-Casher relation 

(gg) = -Vrp(A = 0), (20) 

which relates the quark condensate in the chiral limit to the spectral density of the Dirac 
operator, p(A) = (X^j ^(A — Aj)). For non-zero chemical potential, the eigenvalues are complex 
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and p(A) is a density in the complex plane. How to interpret the Banks-Casher relation in 
that case was explained by Stephanov |^: The quark condensate is related to the resolvent 
of the Dirac operator, [qq) = G{iO) where G{z) = {tT{ip + z)^^) defines a function in the 
complex plane. The vector field G = (ReG, — ImG), can be interpreted as the "electric 
field" generated by the "charge density" p{z). If quarks are condensed this function has a 
discontinuity on the real line. This means that there has to be a finite density of eigenvalues 
on that line. 

If the chemical potential is smaller than the critical value where the baryon density 
becomes non-zero the quark condensate has to be independent of fi. This means that while 
in general eigenvalues will spread in the complex plane, eigenvalues near zero virtuality 
have to stay on the real line. As soon as /i > /ic, some of the small eigenvalues can move 
into the complex plane, but the density of eigenvalues on the real line has to remain finite. 
Only above the critical chemical potential for chiral symmetry restoration can all eigenvalues 
spread out in the complex plane. 

In the remainder of this section we will discuss eigenvalue distributions in the instanton 
liquid model. These distributions were obtained from numerical simulations of the partition 
function (0). For simplicity, we have kept the instanton density (N/V) = lfm~^ fixed. In 
principle, the correct equilibrium density is determined by minimizing the free energy. This 



can be done using the methods discussed in |4^, but is beyond the scope of our exploratory 
study. 

Fig. 1^ shows scatterplots of eigenvalues of the Dirac operator obtained in quenched 
simulations for A^^^ = 2,3 and fi = lfm~^. For = 3, the eigenvalues spread almost 
uniformly over a band in the complex plane. In the case of two colors, there is a finite 
fraction of eigenvalues that stays on the real line. This phenomenon was also observed in 



random matrix models ||51]]. From our discussion above that would seem to indicate that for 
Nc = 2 chiral symmetry remains broken, while for Nc = 3 chiral symmetry is restored. This 
is not the case. Chiral symmetry restoration in the case of three colors is an artefact of the 
quenched approximation. For Nc = 2 the density of eigenvalues on the real line vanishes in 
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the infinite volume limit. We will study this case in more detail in the next section. 

For Nc = 2 the unquenched eigenvalue distribution is qualitatively similar to the 
quenched distribution. This is not true for Nc = 3. In order to calculate the unquenched 
eigenvalue distribution we have binned the eigenvalues and determined the average num- 
ber of eigenvalues per bin from numerical simulations of the partition function The 
unquenched distribution is then obtained by including the phase of the determinant accord- 
ing to ([T9|) . In general, this procedure will give a complex eigenvalue density but in our 
simulations the imaginary part was always consistent with zero. 

In Fig. |]we show the density of eigenvalues for fi = 1.0 fm^^ without the phase (upper 
panel) and with the phase included (lower panel). We observe that if the phase is included, 
the number of eigenvalues away from the real line near A = is reduced. This is clearly the 
correct tendency although with our statistics, we cannot really determine whether there is 
a discontinuity on the real line. We will study the quark condensate more quantitatively in 
the next section. 

VI. CHIRAL SYMMETRY RESTORATION IN THE INSTANTON LIQUID 

In this section we wish to study the quark condensate in unquenched simulations of the 
instanton liquid at finite chemical potential. Let us start with the simpler case = 2, Nf = 
2. In Fig. I we show the quark condensate and the fraction of real eigenvalues as a function 
of the chemical potential for different volumes V = 16, 32, 64 fm^. The instanton density was 
fixed at (N/V) = 1 fm~^ and the quark mass is m = 20 MeV. This mass is somewhat larger 
than the physical quark mass and was chosen to minimize finite volume effects. We observe 
that the quark condensate decreases as /i increases and that chiral symmetry is restored for 
large chemical potentials. While the number of real eigenvalues decreases as the volume is 
increased, the behavior of the quark condensate does not change appreciably. This would 
seem to indicate that for finite quark masses, there is no sharp transition in the infinite 
volume limit. 
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The behavior of the quark condensate and the chiral susceptibihty Xo- = V {{(qqY) — {qq)'^) 
for different quark masses m = 20, 30, 40 MeV is shown in Fig. ^. We have also determined 
the /i = pion mass for these quark masses. We find m^r — 240, 280, 360 MeV, with 
statistical errors AttItt — 20 MeV. The onset fiQ of chiral symmetry restoration is hard to 
determine from the data, but the results appear to be consistent with ~ 171^^/2. We 
also observe that the pseudo-critical chemical potential /ic (determined from the peak of the 
chiral susceptibility) drops with the quark mass. Again the errors are large, but the behavior 
is compatible with /ic ~ i/m- 

For large chemical potential the quark condensate goes to zero but we expect that di- 
quarks are condensed. This idea is based on the fact that SU{2) gauge theory with two 
massless flavors has an extra symmetry that relates the chiral condensate to the 

scalar diquark condensate {e"'^q°'C'')^T2q^) . Here, C is the charge conjugation matrix, a, h are 
color labels and T2 is the anti-symmetric Pauli matrix (the diquark is an iso-singlet). A finite 
quark mass breaks the symmetry and favors quark-anti-quark over diquark condensation. 
For large chemical potential, on the other hand, we expect that diquark condensation is 
favored. This scenario was first discussed in the context of the strong coupling limit by 
Dagotto et al. . However, these authors concluded that chiral symmetry remains broken 
for large /i. In general, this is not correct. For Nj = 2 the diquark condensate is a chiral 
singlet. For more than two flavors, on the other hand, part of the chiral symmetry is broken^. 

Let us now study the more interesting case A^c = 3. We have performed simulations 
with the Dirac operator evaluated in a basis of exact fi ^ zero modes ([l6|) and with the 
schematic Dirac operator (|T7|). The results presented here were obtained by generating 5000 
configurations in a volume V = 16 fm^. In addition to that, some runs were performed with 
V = 32 fm^. To illustrate the procedure, we show trajectories of the quark condensate and 
the phase of the determinant for 400 Monte Carlo sweeps at /i = 0.5 fm~^ in Fig. ^. The 

^Dagotto et al. use staggered fermions, so effectively the number of flavors is a multiple of 4 
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fluctuations are larger as compared to the fi = case, but {e^'^) is still differs significantly 
from zero. If /i > 1 fm~^ the fluctuations are so large that both the real and imaginary parts 
of (e**^) are close to zero. 

Fig. 1^ shows the real parts of the quark condensate and the average phase (e**^) calculated 
with the /i = zero modes. The imaginary part of both quantities is always consistent with 
zero. We display both {qq) and {qq)\\, the quark condensate with and without fluctuations 
of the phase of the determinant taken into account. 

We observe that {qq)\\ drops rapidly with /i. This is the expected artefact that comes 
from neglecting the phase of the determinant. With the phase taken into account, (qq) 
becomes almost independent of the chemical potential for /i < lfm~^. In general we expect 
{qq) to be //-independent for /i < /i^ with /ic — niB/S ^ 1.5 fm^^. In the instanton model 
this is not quite true, because the most tightly bound state that carries baryon number is a 
scalar diquark, so /ic — m/)/2. 

The data show a small but statistically significant rise in {qq) for fi < lfm~^. The result 
is not a finite volume effect because it is also present in the V = 32 fm^ data. A similar 
phenomenon was also observed in random matrix models where it was argued that the 
problem is related to the restriction to a finite number of static modes. 

For chemical potentials larger than /i ~ 1 fm"^ the complex phase becomes very small 
and {qq) has very large error bars. Without the phase, the quark condensate is small and 
fluctuates little. This is due to the fact that the eigenvalues of the Dirac operator are 
pushed farther and farther away from the real line, see Fig. 0d. If the chemical potential 
is very large, /i > 5fm~^, the complex phase factor starts to grow and {qq) can again be 



determined. The reason for this behavior is trivial and was already discussed in Sect. 

If fiMjA is much bigger than the typical overlap matrix element Tj^, the Dirac operator 

becomes anti-hermitean and all eigenvalues are close to the imaginary lineQ. 



^Note that we have performed these simulations without the exponential suppression factor 
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In Fig. ^ we show the results for the Dirac operator evaluated in the basis of exact /i 7^ 
zero modes. Again, {qq)\\ drops rapidly with fi. Including the phase helps to stabilize the 
quark condensate, but there is a definite drop in {qq) even for /i < lfm~^. In this case, the 
effect becomes smaller as we go to larger volumes. We therefore believe that the effect is 
mostly due to the onset of chiral symmetry restoration near ~ lfm~^, smeared out by 
finite volume effects. 

For /i > 1 fm~^ the complex phase factor again becomes very small and we cannot 
perform reliable simulations. However, for fi > 3 iiuT^ the phase factor (e*"^) starts to grow 
and simulations are feasible. In the plasma phase, chemical potentials between (2 — 3) fm~^ 
correspond to a baryon density Ub = 2;U^/ (Stt^) = (0.5 — 1.8) fm~^. In this regime the quark 
condensate is very small and the instanton ensemble is in the chirally restored phase. It is 
interesting to analyze the mechanism that causes the symmetry to be restored. 

For this purpose we show two typical instanton configurations for /i = and /i = 5 fm~^ 
in Fig. 1^. The figures show projections of a four dimensional box into the 0:3 — X4 plane. 
The locations of instantons and anti-instantons are denoted by ± signs, and the strength 
of the fermionic overlap matrix elements Tja_ are indicated by the thickness of the lines 
connecting the instantons. We observe that for /i = instantons are distributed randomly 
and there is no particular pattern to the lines connecting them. This is a sign for chiral 
symmetry breaking. If instantons tend to form large random clusters, then the corresponding 
eigenfunctions will be low-lying delocalized modes. These are the modes that form the quark 
condensate. 

At large chemical potential, instantons tend to line up in the a:4-direction. These chains 
are weakly linked in the spatial direction. The reason why these configurations are favored 

This was done for consistency, because (^) is really a small ^ approximation. Including the 
suppression factor also has the unusual effect that fiMjA never dominates over Tja, because Tja ~ 
p^^ and the factor ex.p(—Nf{p.p)'^) in the instanton rate forces p^^ > p. 
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should be clear from our discussion of the overlap matrix elements in Sect. |TT[ For large fx, 
TjA is big in the temporal direction and oscillates for spacelike separations. There are two 
reasons why these configurations tend to have an almost real determinant. One is the fact 
that for the most attractive orientation, we have U4 = 1 and the overlap matrix element is 
real. The other is that the dominant terms in Tja at large chemical potential are even in /x, 
corresponding to an almost hermitean Dirac operator. 

For the chain-like configurations observed at large chemical potential, chiral symmetry 
is restored because the eigenfunctions tend to be spatially localized. This can be seen from 



Fig. 10 which shows the distribution of participation numbers of the low-lying eigenstates 



for three different simulations atT = /i = 0,T~Tc and /i = [48|, and T = 0, /i > jjc- The 



participation number is a concept borrowed from the theory of localization (Mott- Anderson) 



phase transitions [^. For a normalized eigenvector q, we define the participation number 

p = i/T.h\' (eIc^I' = i)- (21) 

The participation number is a measure of the number of (localized) basis states that con- 
tribute to an eigenstate of the system. We observe that at fi = T = many basis states 
participate in the low- lying eigenstates. In fact, the participation number of these states 
grows with the volume of the system. The finite T transition is associated with the for- 
mation of instanton-anti-instanton molecules and the dominant participation number near 
Tc is two. In the finite /i transition we also observe that the average participation number 
becomes smaller, but the average cluster size is bigger than two. One can check that the 
size of these clusters is dominantly controlled by the inverse temperature (the length of the 
imaginary time direction). We conclude that the system is dominated by oriented instanton 
"polymers" . 

In more physical terms, chiral symmetry restoration is due the presence of valence quarks 
in the system. Valence quarks can saturate the fermionic bonds between instantons. If the 
density of valence quarks is sufficiently large all bonds are saturated, and additional quarks 
can no longer propagate from one instanton to another. As a result, quarks do not acquire a 
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constituent mass and chiral symmetry is restored. This interpretation becomes more obvious 
if one reduces the instanton-induced interaction between quarks to an effective four-fermion 
interaction, and treats this interaction in the mean- field approximation (see, e.g. [0). Chiral 
symmetry breaking is then described in terms of a gap equation for the constituent quark 
mass. At large density chiral symmetry is restored because valence quarks Pauli-block the 
quark loop that drives the formation of a constituent mass. 

VII. SUMMARY 

In summary we have studied instantons and their interactions at finite baryon chemical 
potential. We have emphasized the importance of exact zero modes in the spectrum of the 
Dirac operator even at ^ ^ 0. We have determined matrix elements of the Dirac operator 
between zero modes associated with different instantons. These matrix elements favor the 
propagation of quarks along the positive X4 direction. 

Using these ingredients we have performed a number of exploratory simulations of the 
instanton liquid at finite chemical potential. In the case Nc = 2 and Nf = 2 these simulations 
are straightforward because the theory does not have sign problem. We observe that chiral 
symmetry is restored at fic ~ V^- The high density state is likely to support a diquark 
condensate. 

For Nc = 3 the simulations are complicated by the sign problem. In our simulations we 
focus on the dynamics of (quasi) zero modes. This means that number of degrees of freedom 
in a given volume is much smaller than in a typical lattice simulation and the sign problem 
is not as severe. In practice we can perform simulations in sufficiently large volumes for 
/i < 1 fm~^ and /i > 3 fm~^. For small chemical potential we observe that including the phase 
of the determinant has the correct effect of stabilizing the quark condensate. Qualitatively, 
this is true for both the ensemble generated with the fi = and the correct /i 7^ zero 
modes. Quantitatively, the quark condensate rises for small /i in the first case, while there 
is a tendency towards chiral restoration in the other case. 
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The structure of the chirally restored phase at large ji is completely determined by the 
functional form of the zero mode wave functions. In schematic models of the structure 
of the Dirac operator, the fluctuations in the phase become small because for large fi the 
determinant is dominated by the chemical potential term. Chiral symmetry is restored 
because the large external field pushes all eigenvalues away from zero. In the instanton liquid, 
fluctuations become small because the ensemble is dominated by chain-like configurations for 
which the determinant is almost real. Chiral symmetry is restored because the corresponding 
eigenfunctions are strongly localized. 

Much more work remains to be done in order to understand the nature of the chirally 
restored phase at large chemical potential and the onset of chiral restoration at small In 
particular, we have to study the quark propagator, mesonic correlation functions and the 
thermodynamics of the system at large ii. This should clarify whether there is any tendency 
towards diquark condensation, strange matter formation, etc. 

We should also point out the shortcomings of our procedure. The restriction of the 
determinant to quasi zero modes has enabled us to perform simulations in reasonably large 
volumes but it also introduces certain artefacts. For example, non-zero modes have to 
contribute to the total baryon number of the system. Also, if non-zero modes are neglected 
the quark condensate is not necessarily independent of ii below — ms/'i- Finally, as 
long as there is no real solution to the sign problem, simulations in the interesting regime 

~ will remain impossible. 

Acknowledgements: 1 would like to thank R. Rapp, E. Shuryak, M. Velkovsky and J. 
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APPENDIX A: FERMIONIC OVERLAP MATRIX ELEMENTS 



The two functions /i^2 depend on the instanton-anti-instanton distance {z,Z4), the two 
radii pi^A and the chemical potential /i. If the two radii are equal, it is clear that fi only 
depends on \z\/p, z^/ p and pp. We have checked numerically that if p/ 7^ pa the matrix 
elements depend to a very good accuracy only on the geometric mean p = y/piPA- 

This leaves us with two functions fi{zr, z^, p), where Zr = Zr/p, etc. are dimensionless 
variables. We have determined these functions numerically for many different values of the 
variables and fitted the result to the following parametrizations 
2.0 



fi 



(2.2 + ^2) 



-| (2Z4 + aipzl + a2pz'^ + a^p'^Zij cos{pZr) 



(Al) 



a4/i^ + a^z'^ 



a^zl 



a-jpz4z2j \- asp sm{pzr)j, 



f2 



2.0 



(2 2 + ^2)2 { (^^'' ~ bipzl/Zr + h2^lZ4Zr + hp^ zfj COs{pZr) 

+ (bizl + ^hizlz^ + hpz'^ + bepz'^zl^ sin(pZr) ^ ^ bgp'^z^ } sm 



(A2) 



l{pZr)}, 



where we have dropped the tilde symbol. The functional form of the parametrizations 
was inspired by the free quark propagator (|]). In addition to that we have to respect the 
symmetries fi{z,p) = ^fi{—z,—p), (i=l,2) of the matrix elements and ensure the correct 
behavior for /i — > 0. 

Fitting the parametrization ( [A1| , [X^ ) to the data, we get 



ai 


= 1.334, 


0-2 


= 2.142, 


03 


= 2.872, 




= 0.756, 


as 


= 2.378, 


as 


= 0.893, 


a? 


= 1.060, 


as 


= -0.845, 






hi 


= 2.361, 


b2 


= 1.561, 


bs 


= 0.072, 


W 


= 1.027, 


h 


= 0.818, 


be 


= 0.844, 




= 3.270, 


bs 


= 1.313. 







(A3) 



(A4) 



The quality of the fit can be judged from the results shown in Fig. |I| and appears quite 
satisfactory for our purposes. 
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FIGURES 

FIG. 1. Overlap matrix elements /1.2 as a function of the instanton-antiinstanton separation 
z in units of p. Fig. a) shows /i as a function of 2:4 for Zr = 0, while Fig. b) shows /2 as a function 
of Zr for Z4 = 0. The points are numerical results for the overlap integrals while the lines show 
the parametrization discussed in the appendix. The solid, long-dashed, and short-dashed curves 
correspond to {fip) = 0, 1, 2. 

FIG. 2. Scatterplots of the complex eigenvalues of the Dirac operator from quenched simula- 
tions at /X = 1 fm^^ for Nf. = 2 and N^. = 3. 

FIG. 3. Contourplots of the eigenvalue density in the complex plane from unquenched simu- 
lations at /X = 1 fm^^ without (upper panel) and with (lower panel) the phase of the determinant 
included. The x and y axes show the real and imaginary parts of the eigenvalues. 

FIG. 4. Quark condensate (upper panel) and fraction of real eigenvalues (lower panel) as a 
function of the chemical potential for Nc = 2, Nf = 2 from numerical simulation of the unquenched 
theory. The different curves correspond to different volumes V = 16,32, 64 fm^. 

FIG. 5. Quark condensate (upper panel) and scalar susceptibility (lower panel) as a function 
of the chemical potential in = 2, Nf = 2 QCD. The different curves show the results for different 
quark masses m = 10, 20 MeV with the volume V = 32 fm^ kept fixed. 

FIG. 6. Time history of the quark condensate (both real and imaginary parts) and the cosine 
of the phase of the fermion determinant for an instanton liquid simulation at = 0.5 fm~^. The 
figure shows 400 configurations generated in a volume V = 16 fm'*. 

FIG. 7. Quark condensate (upper panels) and average phase of the fermion determinant (lower 
panels) as a function of the chemical potential for Nc = 3, Nf = 2 from numerical simulation of 
the unquenched theory. The open and solid points show the quark condensate without and with 
the phase included. The calculations shown here were performed in the fi = basis equ. (p!7|). 
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FIG. 8. Same as Fig. with the determinant evaluated in the 7^ basis (|T^). 

FIG. 9. Typical instanton configurations for fi = and fi = 5 fm~^. The plots show projections 
of a four dimensional box into the X3 — X4 plane. The positions of instantons and anti-instantons 
are shown as it symbols. The lines indicate the strength of the fermionic overlap matrix elements 
TiA. 

FIG. 10. Distribution of participation numbers for low-lying eigenvectors in the instanton 
liquid at finite temperature (upper panel), and at finite chemical potential (lower panel). 
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